Dynamical frustration in ANNNI model and annealing 
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Zero temperature quench in the Axial Next Nearest Neighbour Ising (ANNNI) model fails to 
bring it to its ground state for a certain range of values of the frustration parameter re, the ratio of 
the next nearest neighbour antiferromagnetic interaction strength to the nearest neighbour one. We 
apply several annealing methods, both classical and quantum, and observe that the behaviour of the 
residual energy and the order parameter depends on the value of re strongly. Classical or thermal 
annealing is found to be adequate for small values of re. However, neither classical nor quantum 
annealing is effective at values of re close to the fully frustrated point re = 0.5, where the residual 
energy shows a very slow algebraic decay with the number of MCS. 
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I. INTRODUCTION 

Simulated annealing is usually applied to systems with 
frustration, like spin glasses and optimisation problems, 
where the energy landscape is complex with many spu- 
rious minima. There are certain other systems, how- 
ever, which have very simple energy landscape picture 
and ground states, but still the system fails to reach its 
ground state during a energy-lowering dynamical process. 
This situation corresponds to "dynamical frustration" . 
We have specifically considered the case of the axial next 
nearest neighbour (ANNNI) chain, where such a situa- 
tion is encountered. In section II, we elaborate the no- 
tion of dynamic frustration with examples and in section 
III, the dynamics in ANNNI model is discussed in detail. 
The results of application of the classical and quantum 
annealing are discussed in sections IV and V. Summary 
and some concluding comments are given in the last sec- 
tion. 



II. DYNAMIC FRUSTRATION IN ISING 
MODELS 

Quenching dynamics in magnetic systems has been a 
topic of intense research over the last few decades. In 
quenching dynamics, the system has a disordered initial 
configuration corresponding to a high temperature. As 
the temperature is suddenly decreased quite a few inter- 
esting phenomena take place like domain growth [1,2], 
persistence [3-6] etc. 

The Ising model maybe regarded as the simplest model 
describing magnetic properties of many real systems and 
it shows a rich dynamical behaviour with respect to the 
above phenomena. The dynamics of Ising models has 
been extensively studied in lattices of different dimen- 
sions as well as on graphs and networks. 

In dynamical studies, the system is allowed to evolve 
from the initial configuration following a certain prescrip- 
tion and the commonly used dynamical rule at zero tem- 
perature is the Glauber dynamics, i.e., a spin is chosen 



randomly and flipped if it makes the energy lesser, not 
flipped if energy increases and flipped with probability 
1/2 if there is no energy change. 



A A A /ii 



V V V 



A A A A A 



V M 



A A A A 



a A A a 



FIG. 1. Domain coarsening in one dimensional Ising model: 
the spin at the boundary of the left domain wall flips mak- 
ing the two domain walls move closer (Top: earlier picture, 
bottom: later picture). 

The zero temperature deterministic dynamics in Ising 
models can be visualised as the motion of interfaces and 
the domains grow in size as the interfaces annihilate on 
approaching each other (Fig. 1). In one dimension, a zero 
temperature quench of the Ising model ultimately leads 
to the equilibrium configuration, i.e., all spins point up 
(or down). 

In two or higher dimensions, however, the system does 
not always reach equilibrium [7-12]. Such a situation cor- 
responds to dynamical frustration when the system gets 
frozen in a metastable state which does not correspond 
to the ground state. For example, in the two dimen- 
sional lattice (Fig. 2), the dynamics stops at a higher 
energy when the domain walls are straight and appear 
without any corner. The system thus acquires a "striped 
phase" , where the number of stripes is an even number. 
In dimensions higher than two, there maybe other kind 
of frozen states in which the system gets locked. 

This kind of freezing or blocking is also encountered in 
ferromagnetic Ising models on random graphs and small 
world networks where there are finite number of random 
long range bonds. In the random graph, any two spins 
are connected with a finite probability while in the small 
world network, random long range connections occur in 
addition to nearest neighbour links. In these cases, the 
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domain walls may get pinned resulting in a frozen state 
[13-17]. Recently, freezing has been observed on scale- 
free networks also, where, although the system is locked 
in an excited state, the dynamics continues indefinitely 
[18]. 
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FIG. 2. Striped phase in two dimensional Ising model un- 
der zero temperature quenching dynamics - the domain walls 
are straight and no spin flips can occur. 

In the above examples of freezing in Ising models on 
finite dimensional lattices, graphs and networks, a few 
things are to be noted 

(a) The ground state is simple in all the cases 

(b) There is no frustration arising out of the interactions 
in the system. 

In addition, power law scalings with time (e.g., domain 
size ~ i 1 / 2 ) exist for the finite dimensional lattices [19] 
but in the case of networks or random graphs, one has 
an exponential relaxation behaviour consistent with the 
mean field nature of these systems [20-23] . 



III. DYNAMICS IN ANNNI CHAIN 

We will now discuss the dynamics in the axial next 
nearest neighbour Ising (ANNNI) model [24] in which 
there is frustration but no randomness or disorder. The 
ANNNI model in one dimension is described by the 
Hamiltonian 



H 



(1) 



The ground state of this model is well-known: it 
is ferromagnetic for re = J%l J\ < 0.5; antiphase for 
k = J2/J1 > 0.5 and highly frustrated for J%/J\ = 0.5. 
All configurations having domains of size > 2 are ground 
states at the point re = 0.5, which we call the fully frus- 
trated point. 

The dynamics of the ANNNI chain has quite a few 
interesting behaviour [25,26]. When zero temperature 
quenching dynamics dynamics is considered, re = 1.0 
emerges as a dynamical transition point. For re < 1, there 
is no conventional domain coarsening or persistence be- 
haviour. Here, all domains of size 1 immediately vanish 
but domains of size two are stable such that two domain 
walls cannot approach each other and annihilate (Fig. 3). 



Unlike in the Ising model in d = 2, here the domain 
walls are not pinned but can move around keeping their 
number fixed. The energy of the system is constant as 
spin flips occur with zero energy cost. The system thus 
wanders in a subspace of iso-energy metastable states 
forever. As a result, the model also does not show a 
power-law decay in the persistence probability. Interest- 
ingly, there is no special effect of the re = 0.5 point (which 
dictates the static behaviour) on the dynamics. 
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FIG. 3. Dynamics in ANNNI chain for re < 1: All domains 
of size 1 are unstable and they immediately vanish (a). The 
system is left with domains of size > 2 . A domain of size two, 
as shown in (b), is energetically stable as the spins within the 
domain satisfy the antiferromagnetic ordering with their sec- 
ond neighbours and energy contribution from first neighbours 
is zero. As a result, two domain walls have to maintain a mini- 
mum distance and cannot approach each other and annihilate. 

As the domain walls continue moving in the system, 
the number of spin flips at any time becomes a constant 
in time. This constant is independent of the value of re. 
That is expected as this quantity is proportional to the 
number of domain walls. The average number of domain 
walls remaining in the system (per spin) turns out to be 
close to 0.28 [26]. 
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FIG. 4. Constant residual energy of the ANNNI model with 
zero temperature dynamics for re < 1. 

The residual energy E r , defined as the excess energy 
over the ground state energy shows an interesting be- 
haviour with re (Fig. 4). At small values of re the large 
number of domain walls makes the residual energy large. 
As re is increased, E r becomes lesser and at re = 0.5 it is 
zero. This is not surprising however; the configurations 
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with domain sizes > 2 are nothing but the degenerate 
ground states of the k = 0.5 point. The residual energy 
decreases as k = 0.5 is approached from both sides. 

For the sake of completeness, it should be mentioned 
that the ANNNI model shows conventional relaxation be- 
haviour for k > 1 although with a dynamic exponent 
different from that of the nearest neighbour Ising model. 

Since the conventional dynamics does not bring the 
system to the ground state, one has to employ some other 
method to do this. Simulated annealing for the Ising 
model on random graphs has been attempted to "melt" 
the system with success [13]. For small world networks, 
freezing can be got rid off by letting more number of edges 
in the system also [27]. We discuss in the next two sec- 
tion the result of applying different annealing schedules 
to the ANNNI model for various values of n < 1. 



IV. CLASSICAL ANNEALING (CA) 

We have adopted two different schemes for applying 
classical or thermal annealing. 

Scheme A 

This is the conventional scheme where one starts with 
a finite temperature T — To and slowly reduces it accord- 
ing to a linear schedule, such that, at the t-th iteration 
step, 



ground state is also antiphase. The power law remains 
valid for k > 0.5 with an increasing value of a. 



T = T {l-t/r), 



(2) 



where r is the total number of Monte Carlo steps (MCS). 
The final temperature (at t = r) is zero for any starting 
value of T. 

Since for all non-zero temperature, the ANNNI model 
is in a paramagnetic state, one may start with a random 
initial configuration corresponding to To. 

We have calculated the residual energy and the order 
parameter as functions of t; the former is expected to ap- 
proach zero and the latter should increase towards unity 
with larger values of r. 

In Fig. 5, we show the behaviour of E r against r for 
different k values. These simulations have been done for 
a system of 100 spins keeping T = 10. The number of 
configurations n over which averaging has been done de- 
creases with r; starting with n — 1000 for smaller values 
of t, it is decreased to n — 100 for very large r values. 
While for k = 0.2, we find a stretched exponential decay, 
the nature of the curves changes to a power law decay 
for higher values of K. Close to 0.5 it has a very slow 
decay. Corresponding to a power law fit E r ~ t~ q , a 
is very small here, e.g., for k — 0.4, a = 0.08 ± 0.01, 
and for k — 0.6, a — 0.03 ± 0.01. The slowing down of 
the decay pattern of E r seems to depend on the close- 
ness to ft = 0.5 as well as on the nature of the ground 
state. For the various values of k for which the annealing 
scheme has been employed, the slowest decay is observed 
at k — 0.6, when the system is close to K = 0.5 and the 
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FIG. 5. Residual Energy vs. MCS r in the ANNNI model 
for different values of k < 1.0. The decay is slow close to 
k = 0.5. T = 10 here. 

We have also checked the efficiency of this scheme for 
different values of To. Note that the decrease in T is 
made with a slope equal to T /r (eq. 2), so we compare 
the results for different To's by plotting E r against t/Tq. 
For k — 0.4 or 0.6, E r is independent of To (as long as 
To is not very small compared to 1) while for k = 0.2, 
lowering To makes the decay faster at large r. However, 
when To is made smaller than unity the decay does not 
become faster any more. This plot (Fig. 6) also shows 
that for k close to 0.5, the power law behaviour is actu- 
ally valid over a large range of r which is not so apparent 
from Fig. 5. 
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FIG. 6. The residual energy vs the number of MCS (r) for 
different starting values of To. The scaled data shows a good 
collapse for To > 1 for k = 0.4, 0.6. For k = 0.2, the decay is 
fastest for To = 1 (filled circles). The data for k — 0.6 and 
0.2 have been shifted for clarity. 
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Thus E r depends strongly on k: not only does the func- 
tional form change from a stretched exponential to power 
law, there is also a non-universal exponent a which de- 
pends on the value of k. The role of k — 0.5 is felt clearly 
as the annealing is least effective close to this point. 

While estimating the order parameter (OP), it should 
be noted that for k < 0.5, the order parameter is just the 
magnetisation while for n > 0.5, it is the average of the 
four sublattice magnetisations defined as 

m a = S^Zo _1 5 Q+ 4 3 ; a = 1,2,3,4 (3) 

as in [26]. In Fig. 7, we plot the behaviour of the OP with 
t for different values of k. The behaviour of the order 
parameter is also n dependent. The variations with r are 
not smooth: the reason is that decreasing energy is not 
necessarily equivalent to increasing the order parameter. 
However it seems to have a rough power law increase for 
all values of k (Fig. 7). As expected, the growth of the 
order parameter is slowest at n = 0.6. 
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FIG. 7. Order parameter vs. MCS r in the ANNNI model 
for different values of k < 1.0. The growth is slow close to 
k = 0.5. 

Thus we find that this annealing schedule is not very 
effective near n = 0.5 but works well for small values of 
k. It is not possible to detect whether there is a value of 
k for which the behaviour of the residual energy changes 
from stretched exponential to power law from the present 
numerical study. 

Scheme B 

In this scheme, we first let the system evolve from a 
random configuration (corresponding to a high temper- 
ature) using the zero temperature dynamics and then, 
after a few steps, apply scheme A. The difference here 
is, we let the system reach one of the metastable states 
under zero temperature dynamics, heat it to a finite tem- 
perature To which is then gradually decreased. 

In the ANNNI chain, for the first 100 iterations, the 
temperature is kept zero such that when re-heated, the 



system is in one of the metastable iso-energy states. We 
find that this scheme accelerates the decay of the residual 
energy remarkably for small values of k (e.g., k = 0.2) at 
large r. However, for higher values, e.g., k — 0.3, 0.4, 0.6, 
for which scheme A gave a power law decay of E r , the 
results are identical to that of scheme A. Thus near the 
k = 0.5 point, this scheme is also seen to fail to bring the 
system to its static ground state. 

This scheme is highly appropriate for cases where 
the system has a fractional probability to end up in 
a metastable or frozen state not corresponding to the 
ground state. Here it is not possible to predict whether a 
random initial configuration will reach the ground state 
under zero temperature dynamics or not. A good exam- 
ple is the two dimensional Ising model as it reaches a 
frozen state in about 30% cases. It is useless to apply 
scheme A here because 70% of the cases do not require 
any annealing at all. Therefore in order to see whether 
CA is useful, it should be better to apply scheme A to an 
initial configuration of frozen state which can be assumed 
to have evolved from a perfectly random initial state with 
zero temperature Glauber dynamics. Thus effectively it 
has undergone a period of cooling at zero temperature 
and when scheme A is now applied to it, it is equivalent 
to scheme B. 
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FIG. 8. Order parameter (magnetisation) and residual en- 
ergy versus MC time r in the 2-d Ising model under scheme 
B. The value of To is 1 here. The dashed line has slope equal 
to 3.33 x 10 -5 in the log-linear plot. 

We have performed some simulation on a square lattice 
of size L = 40, where the number of stripes s is equal to 
2 or 4. In each case, we find that the behaviour of both 
the residual energy and the magnetisation is exponential 
which means that the scheme works very well in this case. 
In Fig. 8, we show the variation of the residual energy 
and magnetisation for s = 2. The exponential relaxation 
is easily understood, as the thermal perturbation breaks 
the structure of the domain walls and the system is then 
again free to evolve dynamically. The results are iden- 
tical for s — 2 and 4 at large r. It may be mentioned 
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here that for a value of the stripe number s comparable 
to L, the situation is very similar to the ANNNI model 
(as stripe sizes have to be > 2) and then the exponential 
behaviour may no longer be present. However, the prob- 
ability of a large value of s is small and so, we have not 
considered values of s > 4. 



V. QUANTUM ANNEALING (QA) 

Although the classical annealing methods work quite 
well for the ANNNI model for small k, we find that close 
to the k = 0.5 point, it leads to very slow relaxation. 
In several situations, quantum annealing is far more effi- 
cient in decreasing the energy of the system [28,29] and 
we therefore apply this method in the ANNNI model. 
Here instead of thermal fluctuation, quantum fluctuation 
is considered to induce tunnelling to enable the system 
reach the ground state for n < 1. 

The Hamiltonian for the quantum ANNNI chain is : 



H = —JiT<SiSi+i + J2^SiSi+2 — TTiSi 



(4) 



This can be mapped to a 2-dimensional classical model 
[30,31] using the Suzuki- Trotter formula: 

H = —J\HS a ,iS a ,i+l + J2^Sa t iSa t i+2 ~ Jp^S a ,iS a +l,i 

(5) 

where 

J p = -PT/21n(tanh(r/PT)) 

and a denotes the ath row in the Trotter direction. Sub- 
sequently, one can use a linear schedule for the transverse 
field as in [29] and find out the results for E r and the or- 
der parameter. However, it is not possible to make T 
equal to zero in the last step as that would make J p in- 
finite. 

We would first show some curious features of the re- 
sults on applying this method to the ANNNI model and 
then try to justify the results. 

The Suzuki- Trotter mapping is exact for P — > oo but 
it can be a good approximation if PT > 1. One needs to 
find out an optimum value of PT for which E r does not 
change with P. We therefore fix PT and find out E r for 
different values of P following [29]. We first fix PT = 1. 
For small n, e.g., for k = 0.1, the scheme indeed makes 
E r go down with r quite efficiently. However, the value 
of E r for the same r and different values of P shows that 
E r actually increases with P. Thus results for any finite 
P may not be reliable. Even increasing PT to 2, we find 
that this behaviour persists. The reason for this maybe 
that the Suzuki- Trotter mapping works with a non-zero 
temperature for which the system is disordered and is 
always at a high energy state compared to the perfectly 
ordered state and therefore E r does not go to zero for 
large r and P. These results are shown in Fig. 9. 
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FIG. 9. The quantum annealing scheme does not appear to 
work well for the ANNNI chain. For example at k = 0.1, an 
optimum value of PT is difficult to find out as E r increases 
with P where PT is fixed at 1 or 2 at any value of r. 
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FIG. 10. The variation of residual energy and order param- 
eter for k = 0.4 against r are shown for 3 different values of 
P. E r actually remains constant with r and P while the OP 
(m(t)) seems to fluctuate around a constant value. 



In case of a value of k close to 0.5, e.g., k = 0.4, we find 
that E r remains virtually a constant for all P values when 
PT = 1 which apparently implies that PT = 1 is a good 
optimum value. However, E r actually remains a constant 
for all t values as well showing that it does not relax at 
all. Thus here too the quantum annealing method will 
not work well. The reason is again because a non-zero 
temperature of the system has been used. However, the 
manifestation of this non-zero temperature is different for 
small and higher values of n. 

On hindsight, it may appear that quantum annealing 
is a redundant exercise in this case. However, it is inter- 
esting to find out how the redundancy makes itself known 
for different values of k in different ways. 
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VI. SUMMARY AND CONCLUSIONS 

In summary, we have shown that in systems with dy- 
namic frustration, simulated annealing can be applied 
which gives results according to the nature of the sys- 
tem. For the ANNNI model, which has a competing in- 
teraction leading to frustration, (but well defined ground 
states with trivial degeneracy for k 0.5) classical an- 
nealing seems to work well for values of k close to k = 0. 

We have applied a different scheme of thermal anneal- 
ing where the system is heated after an initial period of 
cooling. The results remain same in case of the ANNNI 
model when the frustration parameter is appreciable. 
But this method is useful in case of some other models 
where the conventional scheme is not very handy, e.g., the 
two dimensional Ising model. The better effectiveness of 
the annealing scheme for the two dimensional model in 
comparison to the ANNNI model may be attributed to 
the fact that there is no frustration in the former. The 
frustration present in the ANNNI model especially near 
k = 0.5 makes it unresponsive to the thermal anneal- 
ing while the unfrustrated two dimensional Ising model 
responds fast to the thermal fluctuation. 

In case of the ANNNI model, it is also interesting to 
note that while in the dynamical studies it was shown 
that the k = 0.5 point hardly has any role to play, and 
dynamical quantities like persistent probability, number 
of spins flipped at any time etc. were k independent, 
things become strongly k dependent under any anneal- 
ing schedule with non-universal exponents governing the 
power law decays. The effect of k — 0.5 is also felt as the 
annealing fails to make any impact close to it. 

The application of a quantum annealing method with 
non-zero temperature also fails to "melt" the system to 
its ground states and its effect is differently manifested 
for small and large values of ks. 

In fact we find that much is left to be done for a suc- 
cessful annealing programme near k = 0.5, e.g., one may 
attempt a zero-temperature quantum annealing schedule. 

Even for the classical annealing case, the change in be- 
haviour of the residual energy as k is varied requires to 
be studied more intricately and possibly for larger syatem 
sizes. Our present study may act as a guideline for such 
future research work. 

Lastly, we note that the variations of the residual en- 
ergy does not follow a Huse-Fisher [32] type scaling in 
any parameter range for the ANNNI model. 
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